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We demonstrate that a highly frustrated anisotropic Josephson junction array(JJA) on a square 

QQ , lattice exhibits a zero-temperature jamming transition, which shares much in common with those in 

granular systems. Anisotropy of the Josephson couplings along the horizontal and vertical directions 

plays roles similar to normal load or density in granular systems. We studied numerically static and 

£N| , dynamic response of the system against shear, i. e. injection of external electric current at zero 

temperature. Current-voltage curves at various strength of the anisotropy exhibit universal scaling 

features around the jamming point much as do the flow curves in granular rheology, shear-stress 

vs shear-rate. It turns out that at zero temperature the jamming transition occurs right at the 

^H . isotropic coupling and anisotropic JJA behaves as an exotic fragile vortex matter, it behaves as 

^s ' superconductor (vortex glass) into one direction while normal conductor (vortex liquid) into the 

^vj t other direction even at zero temperature. Furthermore we find a variant of the theoretical model for 

the anisotropic JJA quantitatively reproduces universal master flow-curves of the granular systems. 

i I 1 Our results suggest an unexpected common paradigm stretching over seemingly unrelated fields - 

i-C ' the rheology of soft materials and superconductivity. 

<a : 

Physics continues to thrive on analogy [l[ . Rheological properties of matters [2[ and electric transport properties 

of superconductors [3[ exhibit intriguing analogies. The flow curves in rheology, the shear-stress vs the shear-rate, 

& ■ correspond to the current- voltage curves in superconductors |J| . Exploring this analogy further we here demonstrate 

j/2 ' by computer simulations that there exist the zero-temperature jam ming transitions and glassy non-linear rheology, 

+J . originally found in granular and other materials [a, |g, 0, M, H LLiA HH LL2l [Dl, LLJ, UM, UM i m a class of highly frustrated 
anisotropic Josephson junctions arrays (JJA) on a square lattice. Our key observation is that anisotropy of Josephson 
coupling plays the role of normal load or density in granular system such that a jamming transition takes place in the 
limit of isotropic Josephson coupling at zero temperature. Combined with accumulating evidences that the (vortex) 
liquid-glass transition occurs at zero temperature in isotropic JJA, our result provides a strong evidence that the 
isotropic coupling point at zero temperature is an ideal example of the so called J-point (Jamming point) [fj( and 
fi | that the anisotropic JJA is a promising system which allows explorations of both (athermal) unjamming- jamming 
transition and (thermal) liquid-glass transitions in a unified manner in a single system as originally proposed in 
the context of granular and glassy materials [6|. Furthermore, we show that a variant of the original JJA model 
emphasizing elastic nature in a particular direction can quantitatively reproduce scaling features of granular jamming 

pr\ ' transitions observed near the J-point. 

■<^j- | The JJA [3|, [jj], 13 i s a network of superconducting islands as depicted in Fig. 1 (a). The phase of superconducting 

^sO • order parameter of the islands at site i, 9i, is coupled to its nearest neighbors by the Josephson junctions. External 
transverse magnetic field B thread the cells in the forms of flux lines each of which carrying a flux quanta </>o . On 

r~^ . average each unit cell of the square lattice carries / = Ba 2 /4>o flux lines. A flux line threading a unit cell induces a 
vortex of the phases around the cell much as a dislocation in a crystal. Here Aij is the vector potential. The static 
and dynamic properties of the JJA under transverse magnetic field are known to be described to a good accuracy by 
the energy term associated with the Josephson couplings @, [l7| , 
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H = J2 K H E i3 (h - 9j - Mj ) E l3 (VO = - cos(V>) (1) 

<i,3> 

where Ky is the strength of the Josephson coupling. 

In mid 80's a tantalizing possibility of a superconducting glass in the JJA has been raised by Halsey [19(: if the 
transverse magnetic field is tuned in such a way that the number density / of the vorticies takes irrational values, a 
glassy state may be realized at low temperatures because the vorticies may not be able to form periodic structures 
called vortex lattices which are analogous to ordered structures of dislocations in Frank-Kasper phases. It has been 
argued that the frustration due to the gauge field like Ay in Eq. ([T]) mimics geometric frustration in structural glasses 
(see |2fj| | for a review on the perspective on the frustration-based point of view on the glass transition). Indeed 
equilibrium relaxations were similar to the primary relaxation observed in typical fragile supercooled liquids [2li |22| . 
Now recent studies appear to convincingly suggest that the putative (the rmal) liquid-glass transition is actually taking 
place only at zero temperature exhibiting diverging length scale(s) j23l [24l. [25f. We note that all these observations 
are made on systems in which the Josephson couplings K^ are isotropic. 
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FIG. 1: Highly frustrated anisotropic JJA and granular materials. We consider an array of Josephson junctions on a 
square lattice with linear size L (in the unit of lattice spacing a) (a). The ladder version of the JJA model (c) can be mapped 
onto a friction model (d). The control parameter corresponding to the number density of granular particles <f> or normal load 
(f) is the strength of the anisotropy of the Josephson coupling A (e). 



As stated at the beginning our key observation is that the anisotropy of the Josephson coupling K ij is a relevant 
variable and plays the role similar to the particle normal load or density in the jamming of granular materials [7J 
(Fig. 1 (a)(b)). Wc parametrize the anisotropic coupling as, 



A", 



1 for (i,j) directed along the x axis 
A for (i,j) directed along the y axis 



(2) 



Such an anisotropic JJA can be realized experimentally by controlling the width of the junctions [26|]. Let us also note 
that similar anisotropic couplings arise naturally also in cuprate high-Tc superconductors 27) and charge-density-wave 
systems [2a| . 

Another key observation is that the model Eq. ([I]) can be slightly modified to build an effective Eulerian model 
for rheology of granular systems under horizontal shear as shown in Fig. 1 (b). We assume that particles move 
predominantly parallel to the direction of driving (a; axis) confined within the horizontal layers without passing over 
others in the same layer. Then the variable 9i can be interpreted as a phase variable [J| by which the number density 
of particles within the i-th cell is described, for instance, as pi = (1/2)(1 + cos(#j)). Here the size a of the cell 
corresponds to the typical scale of particles and their mutual distances [4y|. The sinusoidal intra-layer couplings 
Eij(ip) = — cos(-0) in Eq. {1} are replaced by elastic couplings Eij(ip) = ip 2 /2 while the sinusoidal form is kept for the 
inter-layer couplings to allow phase slips between different layers. Let us call such a model as semi-elastic model. Here 
the vorticies represent dislocations. We assume that t he g eometrical frustration induced by the gauge field Aij mimics 
real frustrations in granular and other glassy systems [201 ] . In the absence of the frustration (A^ = 0) the semi-clastic 



model exhibits a non-linear rhcology associated with a Kostcrlitz-Thouless transition at finite temperatures [4J. 
The dynamics of the models can be described by the equation of motion, 

d6i dv,i dH 

-dF = v *> m n = -w l +F * (3) 

Here the frictional force Fi is given by Fi = — 7^,-("i~ v j) with the summation taken over the 4 nearest neighbours of 
i. This equation of motion is nothing but the standard resistively and capacitively shunted junction (RCSJ) dynamics 
[3, [2J] , which can also be viewed as a model for rheology of the layered systems [4[ ■ For simplicity we choose mass 
(capacitance) m = 1, damping constant (resistance) 7 = 1. Here we focus only on the zero temperature dynamics. 
Appropriate thermal noise can be added to Eq. ([3]) for finite temperature dynamics. For the semi-elastic model we 
assumed two different types of constitutive relations for the frictional force Fi in Eq. $5§: (1) Newtonian viscous 
friction Fi = — ^2/Avi — Vj) and (2) Bagnold's friction [3a . |39| ] Fi — — J^. \vi — Vj\(vi — Vj) where the sum is taken 
over the nearest neighbours on the two adjacent layers. 

The anisotropic coupling Eq. ^ can be motivated by recalling a well known problem in the science of friction. A 
class of friction models related to the Frenkel-Kontorova model [29| is known to exhibit the so called Aubry's transition 
[30( which is a kind of jam ming transition at zero temperature. In Fig. 1 (d) we display a friction model proposed by 
Matsukawa and Fukuyama [31| which consists of two layers of atoms representing surfaces of two different solids. In 
general the ratio (winding number) / = l/a of the mean atomic spacings on the two different materials takes irrational 
values. The atoms in the same layers are connected to each other by springs while those on different layers interact 
with each other via short-ranged interactions of strength A, which mimic the normal load. In the weak coupling 
regime A < A c , the two chains of atoms slide smoothly with respect to each other thanks to the incommensurability. 
On the other hand in the strong coupling regime A > A c the system is pinned into amor pho us metastable states and a 
finite static frictional force or yield stress emerges. Note that an anisotropic JJ ladder [32| shown in Fig. 1 (c) which 
consist of two horizontal layers can be viewed as an Eulerian formulation of the friction model. The irrational winding 
number / can be identified with the irrational number density / of vorticies in a unit cell in the JJ ladder. 

An important consequence of the anisotropic coupling Eq. ([2]) in 2 (and higher) dimensions is that the effective 
repulsive long-ranged interactions between the vorticies become anisotropic. For A < 1 the vorticies will tend to align 
vertically since the repulsive force is stronger along the x axis, which make it much harder for the vorticies to move 
along the y axis, i.e. the direction with weaker coupling. Of course the situation becomes reversed for A > 1. 

The rigidity of the system can be probed by applying an external current just as external shear stress is applied on a 
solid. What corresponds to the shear stress T, xy along the x axis (Fig. 1 (b)) is the vertical external electric current J y 
(Fig. 1 (a)) |J|. Then vorticies (dislocations) are driven along x axis by the Lorentz force. The resultant electric field 
E y which is proportional to the average velocity of the vorticies corresponds to the shear-rate T xy which measures 
the rate of plastic deformations in rheology. If the vorticies don't move significantly resisting against the Lorentz 
force, the energy dissipation is negligible and the system remains macroscopically superconducting. In practice, we 
apply shear to the system by forcing the top and bottom layers (walls) to move along the opposite directions at 
constant velocities. We measure the resultant electric field E y (shear-rate T xy ) defined as the slope of phase velocity 
Vi developed in the system along the y axis. Electric current flowing through a junction from site i to j is defined 
as sin(^i — 9j — Aij). The currents running through the junctions parallel to x and y axes correspond to the shear 
stresses T, xx and Y* xy respectively in rhcology. 

By construction of the system, static and dynamic resposcs to J x at anisotropy A at T — is just the same as those 
to J y with A' = 1/A. Thus in the following we only display results of resposcs to J y . 

We numerically solved the equation of motion Eq. (|3|) by the 4th order Rungc-Kutta method [47| . Periodic boundary 
condition is imposed along the x axis only. For a given irrational vortex density / we used its rational approximations 
p/q with integer p and q in systems of sizes L = nq (with n = 1 or 2). To explore larger length/time scales we use 
systematically better approximants to prevent commensurability (or matching) effects. (See APPENDIX [Aj Before 
starting measurements, we checked that the velocity profile becomes linear into the y axis without shear-bands and 
that observed quantities do not depend on the prior shear histories. 

Shown in Fig. 2 are snapshots of the vorticies and local currents under shear of the anisotropic JJA. Chains of 
electric currents reminiscent of "force chains" [5|, |7| in granular materials can be noticed. The confi gur ations at the 
isotropic point A = 1 appear to manifest the diagonal stripe structures found in the ground states [331 I34J [35[ . As 
expected the vorticies and the electric currents flowing along the trains of vorticies tend to align into the direction 
with weaker coupling. This observation strongly suggests that a jamming transition takes place at the isotropic point 
A c = 1. The system behaves as a fluid (unjammed phase) for A < A c and amorphous solid (jammed phase) for A > A c 
with respect to J y as depicted in Fig. 1 (c). Furthermore it is interesting to note that the jammed state is inevitably 
fragile in somewhat similar sense as proposed in the context of granular matters [5| : the system with a given A can 
resist against shear only into one direction (i. e. superconducting). Thus we may call such a state of matter as fragile 
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FIG. 2: Snapshots of vorticies and local electric-currents in the steady states under shear at zero temperature. 

Vorticies (a-c) tend to align into the direction with weaker coupling. Local electric currents (d-i) exhibit chain-like configurations 
which are reminiscent of the "force chains" observed in granular materials. The "current-chains" tend to percolate into the 
direction with weaker coupling. The system size is L = 42 with / = 8/21 which approximates (3 — v / 5)/2 = 0.3819... The 
system is sheared such that the electric field E v (shear rate) is 0.006. 



vortex matter in the same spirit of [5j. The qualitative features are essentially the same in the semi-elastic model 
except that vorticies move only into x axis in the latter model. 

The current-voltage curves obtained at different values of the coupling A are displayed in Fig. 3 (a). At stronger 
coupling A > 1 it appears that a non-zero critical current J C (X) = liniE^o J(E, A) exists, which becomes larger with 
increasing A. This means that the Lorentz force does not drive the vorticies significantly so that the system remains 
macroscopically superconducting along the y axis at strong enough coupling A. The finite critical current corresponds 
to the yield stress E c in rheology. The disordered configurations of vorticies shown in Fig. 2 suggests that the system 
is an amorphous glassy state of vorticies. On the other hand, at smaller coupling A < 1 and low enough E the Ohm's 
law J = <j{X)E holds with finite linear conductivity cr(A) which becomes larger with increasing A. Thus the vorticies 
can flow easily producing significant energy dissipation at weak enough coupling A. At the isotropic point A = 1, we 
find a power law J oc E l ~ a with 1 — a = 0.34(3). This corresponds to the so called shear-thinning behaviour (a > 0) 
in rheology 0]. 

The above results strongly indicate that A c = 1 is the critical point of a 2nd order phase transition at zero 
temperature. This is supported by a good scaling collapse of the data onto a master curve as shown in Fig. 3 (d). 
Our scaling ansatz is similar in spirit to the ones used for the usual normal-to-superconducting phase transition at 
finite temperatures [3(| [37| which can also be reinterpreted in the context of rheology [4[ , 
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The scaling function (master flow curve) is expected to behave asymptotically as J(x) ex x for small enough x in the 
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FIG. 3: Current-voltage curves and flow curves. The left panels (a,b,c) show the current- voltage of the JJA (a), flow 
curves of the semi-elastic model with Newtonian friction (b) and that with Bagnold friction (c). The right panels (d,e,f) 
are corresponding scaling plots. The master-curves obtained in numerical simulations of 2 dimensional granular systems [4J| 
with coefficient of restitution equal to/smaller than 1 are included in (e) and (f) respectively. In the right panels (a,b,c) the 
strength of the inter-layer coupling A is varied as A = 0.625, 0.75, 0.875, 0.95, 1.0, 1.050, 1.125 and 1.5 from the bottom to the 
top curves. The dotted lines are power law fits to the curve at A c = 1 which yields 1 — a = /3/A = 0.34(3) (JJA), 0.43(6) 
(semi-elastic-Newtonian) and 0.65(4) (semi-elastic-Bagnold) . The scaling plots in the right panels (d,e,f) show J oc J/ (A— A c )^ 

vs E oc E/(\ — Ac) and E oc E/(A — \ c ) vs F oc T/(A — A c ) with A c = 1. The resultant values of the exponent A which 
give best scaling collapse are 3.5 (JJA), 4. 2 (semi-elastic-Newtonian) and 2.4 (semi-elastic-Bagnold). In the left panels (a,b,c) 
/ = 34/89 is used which approximates (3 — v / 5)/2. In the right panels (d,e,f) data of / which approximate 2 — \/3 (d,e) and 
(5 — \/T7)/2 (f) are also included. We have checked that finite size effects and commensurability effects are not significant within 
the range of shear-rates used here. In the semi-elastic-Bagnold model (f ) the system size is limited to avoid strong shear-banding 



effects. The scaling functions of the granular systems (e,f) are obtained by plotting E oc E/(<^> — <f> c )P vs f oc f /{4> 
with (f> c — 0.8415 (random close packing density in two dimensions) and appropriate exponents /3' and A'. 



Ohmic phase (A < A c ) and lim^^o J(%) —* const in the superconducting phase (A > A c ). The scaling ansatz Eq. (j4|) 
implies 1) the linear conductivity diverges as cr(A) oc (A — X c )~' l - A ~^ for A — > A~, 2) the critical current vanishes as 
J C (A) ex (A — A c )^ for A — > A+ and 3) the critical behaviour J(x) oc x x ~ a = x 13 ^ sets-in for large x. Here we used 
the notations reflecting the analogy with the equilibrium critical behaviour of ferro-magnets under magnetic field as 
noticed by Wolf, Gubser and Imry [361 ] : the shear plays the role of symmetry breaking field like the magnetic field 
and the critical current emerges as an order parameter like the magnetization (see [ll| for a similar argument in the 
context of rheology). As shown in Fig. 3, we find our scaling ansatz works well with A c = 1. We have checked that 
the universality does not depend on the use of different irrational values of / as demonstrated in Fig. 3 (d). 
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FIG. 4: Static helicity (shear) modulus in the anisotropic JJA. The helicity modulus Y y is obtained as equilibrium 
ensemble averages of the susceptibilities against twists along the y axis, (a) Temperature T and size L dependence of the 
helicity modulus Y y and (b) its extrapolation to T — > 0. In (b) helicity modulus at T — is also shown. The values of / used 
here approximate (3 — V&)/2. The periodic boundary condition is imposed along both x and y aixes on the system of size 
N — L x L. The standard simulated annealing method is used to generate the equilibrium ensemble (appendix [B}. 



Let us note that the critial current discussed above J C (A) = lim.£_,o J{E, A) corresponds to the dynamical yield 
stress in rhcology. On the other hand one can also define the quasi-static critical current (s) needed to move out of 
a generic metastable state (or the ground state), corresponding to the static yield stress in rhcology, as studied by 
Teitel and Jayaprakash [13] on the JJA. In practice one can consider athermal quasi-static processes similar to those 
used in some recent studies on amorphous solids [40, |4l|, [42| : starting from a metastable state reached from a random 
initial configuration by a quncch, i. c. deterministic energy descent process, the system is subjected to externally 
induced uniform strain 6 (See Eq. (|B1[) ) which is increased step by step. The system relaxes down to an energy 
minimum by the energy descent process after each small increment of S. As the result one finds that the current J (6) 
is a sawtooth-like function of 5: piecewise linear lines corresponding to elastic deformations broken by yield points 
at plastic events, as observed in amorphous solids [4fJ, |4j], |42j . A critical current can be defined as the value of the 
current just before reaching a yield point. In [lj| Halsey has found that typical value of such quasi-static critical 
current J ( ! tatlc (l) is finite in the isotropic system A = 1. Moreover we found that J^ tatlc (A) varies smoothly with A and 
remains finite even in the unjammed phase A < 1 found above [43[. Our data of J(E, A) in the flow curves becomes 
smaller than J^ tatlc (A) meaning that the quasi-static current needed to move out of a generic metastable state and 
dymamic critical currents are distinct in the present system. Quite interestingly we observed that the flow curves 
become strongly dependent on strain histories below J^ tatlc (A). In practice we had to use an annealing procedure to 
obtain stationary data: decrease E (shear rate) very slowly down to the target one. Slower shear rates are needed to 
investigate smaller current J (shear stress) regions. Note that no annealing is performed in the athermal quasi-static 
process discussed above. These observations suggest ruggedness of the energy landscape of the present system. 

We also investigated static response to shear (Sec APPENDIX [B] for the details.). As shown in Fig. 4, the static 
helicity (shear) modulus is very sensitive to the anisotropy. The figure shows that the helicity modulus remains zero 
down to T — > for A < 1 and becomes finite for A > 1. Thus at T = the isotropic point appears to be the 
critical point A c = 1 being consistent with the dynamic response discussed so far. Here let us recall again that just 
by symmetry, helicity modulus Y x /\ at anisotropy A is identical to Y y /X' at A' = 1/A. Then the fact that A c = 1 
means that this system is quite exotic: fragile vortex matter which behaves as a solid with respect to shear along one 
direction but liquid for the other direction. 

Another remarkable feature is the difference of the static helicity modulus between T = and T — > limit which 
can be seen in Fig. 4 b). The helicity modulus at T — only reflects local stability of energy minima. This difference 
suggests the existence of certain softmodes with vanishingly small energy gaps as discussed in APPENDIX [Bj This 
observation suggests that the vortex liquid behaviour at T = is realized by some non-trivial softmodes due to 
frustrations. 

Based on the above results we obtain the jamming phase diagram of the anisotropic JJA as shown in Fig. 1 (c) 
which is surprisingly similar to that of granular systems shown in Fig. 1 (f ) [g, 0, @] • 

Now let us turn to the semi-elastic model. In Fig. 3 (b,c), we display the flow curves of the semi-elastic model. 
The shear-stress due to the inter-layer coupling terms also obeys the Newtonian or Bagnold scaling for small A at low 
enough I\ Flow curves at large A suggests existence of non-zero yield stresses S C (A) = linip^ S(r, A). We find again 



a power law behaviour E oc (T) 1 a with a > (shear-thinning) at A = 1. Indeed the scaling ansatz Eq. (|4]) (with 
E — > r and J — > S) works well again assuming A c = 1 and the universality does not depend on the different irrational 
values of / (Fig. 3 (e,f)). 

Recent numerical simulations of granular materials with/ without strong dissipation at the particle level (coefficient 
of restitution smaller than/equal to 1) have found Bagnold/Ncwtonian scalings in the fluid phase and different critical 
exponents pjj, LL2|, LUa, LLJ, H]- Quite interestingly the values of the shear-thinning exponent 1 — a found in our 
semi-elastic model with Bagnold/Newtonian frictions are 0.63 and 0.42 respectively in agreement with the exponents 
of the corresponding two-dimensional granular systems with repulsive linear spring forces between the particles. For 
a comparison master flow curves of the granular systems |44j ar e displayed in Fig. 3 (e,f). Quite remarkably the 
functional forms of the master flow-curves themselves agree very well. 

In the present paper we focused on the responses of the systems to shear at zero temperature. Recent studies at 
the isotropic point A = 1 at finite temperatures suggest critical behaviour with T —* 0, i. e. T C (X = 1) = with 
diverging length scales [23l . |24| , |25| | . If this would be confirmed, our system would provide a fascinating example where 
both the (thermal) liquid-glass transition and the (athermal) unjamming-jamming transition take place at the same 
thermodynamic point, T = and A c = 1, demonstrating deep connection between the two transitions. At least in this 
system, the jamming and glass transitions appears to be the two sides of a coin. It would be interesting to further 
explore this connection to make this statement more substantial. A related interesting question is whether or not the 
jamming (glass) phase survives at finite temperatures. Then the possibilities arc 1) the jamming point at A = 1 at 
T = is an isolated critical point or 2) a critical line T C (A) starts from the jamming point as shown in Fig. 1 (e). Our 
preliminary study points to the latter possibility. 

We emphasize that the jamming transition here is purely due to geometrical frustration which is free from any 
quenched disorder in sharp contrast to the conventional vortex glasses [37| for which presence of random pinning 
centers are crucial. This is a much awaited, concrete example of a jamming-glass transition purely due to geometrical 
frustration [2fJ]. It is tempting to speculate that similar phenomena may exist in frustrated magnets such as antifcr- 
romagnets on triangular, kagome and pyrochlore lattices [45| . We also note that it will be important and interesting 
to clarify how quenched disorders, which may not be completely avoided in experimental JJAs, come into play. 
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APPENDIX A: RATIONAL APPROXIMANTS AND COMMENSURABILITY EFFECTS 

The purpose of the present paper is to analyzed physical properties of the anisotropic JJA with irrrational vortex 
density /. However in order to use the periodic boundary condition we had to use rational approximants for a given 
irrational number. Here we explain how commensurability or matching effects emerge and explain how we avoid them 
in the present study. 

For convenience we considered a family of irrational numbers called quadratic irrationals. A well known example 
is the golden mean (1 + v5)/2. Their rational approximants /„ can be generated systematically by solving a simple 
recursion formulae, 

fn = -^-— Aa n+1 = Ba n + Ca„_i (Al) 

where A, B and C are integer coefficients and ao, a\ are certain integers. It is easy to verify that the f n s converge to 
an irrational number / = lim rwoo /„ which is a solution of a quadratic equation C/ 2 + Bf — A = 0. For instance the 
golden mean (1 + V5)/2 can be obtained using the Fibonacci numbers a n which satisfy the above recursion formula 
with A = B = C = 1 and 02 = ai = 1. Examples of rational approximants /„ are shown in Fig. 5 (a) which 
approximate (3 — V5)/2. They are generated by solving Eq. (|A1[) recursively with A = 1, B = 3, C = — 1 and 
02=01 = 1. It can be seen that the approximation becomes better such that \f — f n \ becomes smaller as n increases. 
In Fig. 5 (b) flow curves of the anisotropic JJA model at T = with the rational vortex densities /„ arc shown. 
Apparently flow curves of a given anisotropy A converge to a limiting curve as n is increased. We regard the latter 
as the flow curve of / = (3 — v5)/2 at anisotropy A. It can also be seen that flow curves of a given approximant f n 
closely follow the limiting curve at large enough electric field E (shear rate F) and deviate from it at lower E. In 
the n dpendent branches the current J(E) (shear force S) tends to saturate to some finite values in E — > limit, 
i. e. critical current J C (A, f n ) (yield stress £ c ) which decreases with increasing n being consistent with the prediction 



by Teitel and Jayaprakash [Tj] (see also [23]). The latter behaviour suggests that a periodic vortex lattice (crystal) 
J17I , Ha , [33l . [34| associated with a given rational vortex density /„ is formed and that the latter dictates physical 
properties of the system at length/time scales larger than its lattice spacing. Thus we expect the so called 'Bingham 
fluid' behaviour (fluid with finite yield stress) cannot be avoided for any small A for fixed n and that the genuine 
fluid phase at T = is realized only for truly irrational vortex densities. On the other hand the above results shown 
in Fig. 5 (b) suggest that physical properties at short enough length/time scales, which are n independent, reflect 
those properties of irrational /. Thus our strategy in the present work is to choose large enough n such that the n 
dependency do not become relvant within the range of E (shear rate T) we choose to work on. 



(a) 



0.45 



0.4 



0.35 



0.3 



(b) io° 



J 



X=1.33 



3 4 5 6 7 

n 



\=1.Q 



E 



^.=0.75 

1 = 1/3 ■ 

2/5 

3/8 
5/13 
8/21 



FIG. 5: Rational approximants and commensurability effects on the flow curves. (a) Approximates f„ for 

/ = (3 - v / 5)/2. (b) Data of the flow curves of the anisotropic JJA at T = with / = 1/3 (L = 42), 2/5 (40), 3/8(40), 5/13 
(39) and 8/21 (42) are shown. The dotted straight line represents the power law J oc E 035 at the critical point. 



APPENDIX B: STATIC RESPONSE WITH RESPECT TO SHEAR 



To study static response to shear along, say y axis, it is useful to consider a modified Hamiltonian, 



*=- £ 



cos 



<*,j>l|ex 



,)-A Y^ cos(<t>ij-5), 

<i,J>||e v 



*>ij = 9i — Oj — M} 



(Bl) 



with periodic boundary conditions along both x and y axes. Here e x and e y are unit vectors along x and y axes. 
This is equivalent to consider a system with twisted boundary condition (shear) with total phase difference Lb forced 
across the system along the y axis [ij} ■ 

In equilbrium at temperature T the free-energy F(6) of the system under shear strain S can be defined. Then the 
static helicity modulus Y y is defined as, 



Y„ 



d 2 F 



<i,j>\\e y 

~ X/3 Yl Yl (( sin 0fe) sin (0w)) - (sin(<?3ij)}(sin(0feO}) 

<i,j>||e B <fc,;>||e„ 



(B2) 



where j3 = 1/ (ksT) and (. . .) are equilibrium thermal averages at finite temperature T. The 1st term on the r. h. s 
reflects direct elastic response around energy minima with respect to a small externally induced shear strain. On the 
other hand, the 2nd term reflects relaxation of the system against the external strain at finite temperatures. 

Here the distinction between T = and T — > is important. At T = only the 1st term exists. However, the 
contribution of the 2nd term can remain, in principle, in T — ► limit if the strength of the thermal fluctuation of the 
current sin(0jj) is 0(T). Such a situation can arise if there are soft modes with vanishingly small energy gap such 
that they remain thermally active at arbitrarily low temperatures. 
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To obtain the equilibrium ensemble to evaluate the hclicity modulus at 5 — 0, we performed simulations of relax- 
ational dynamics by numerically solving the Langevin equation 

with the Hamiltonian H given in Eq. {T]) and &(£) being Gaussian noise with zero mean and unit variance, by the 
2nd order Runge-Kutta method. The system is cooled with cooling rates dT/dt = 10~ 10 — 10™ 9 starting from initial 
temperature at T = 0.3 — 1.0. We have checked that the cooling rate is slow enough by comparing with the results of 
4 times faster cooling rate. 



